StupidButtons
WangJY
2022-11-08
clean data
library(dplyr)
library(readxl)
library(tidyverse)
# nie<-read_xlsx("Data-NIE and WANG(1).xlsx",sheet = 1)
# pra<-read_xlsx("Mobile buttons.xlsx",sheet = 1)
# lbw<-read_xlsx("Data Lbw.xlsx",sheet = 1)
# alldata<-bind_rows(nie,pra,lbw)
# write.csv(alldata,row.names=F,"alldata.csv")
alldata<-read_csv("alldata.csv")
# alldata[alldata$Species=="Rhododendron leptothrium","Species"] <- "Rhododendron_leptothrium"
# write.csv(alldata,"alldata.csv",row.names=F)
summary(alldata)## ForestType No.Plot Elevation Slope
## Length:104 Min. :1.00 Min. :500 Min. :0.500
## Class :character 1st Qu.:1.00 1st Qu.:500 1st Qu.:1.600
## Mode :character Median :2.00 Median :600 Median :4.000
## Mean :2.01 Mean :675 Mean :3.945
## 3rd Qu.:3.00 3rd Qu.:800 3rd Qu.:6.000
## Max. :3.00 Max. :900 Max. :8.200
## DistRiver CoordinateX CoordinateY Species
## Min. : 1.00 Min. : 1.000 Min. :0.0500 Length:104
## 1st Qu.: 6.80 1st Qu.: 2.250 1st Qu.:0.1000 Class :character
## Median :11.00 Median : 4.000 Median :0.2125 Mode :character
## Mean :10.35 Mean : 4.894 Mean :0.2719
## 3rd Qu.:13.00 3rd Qu.: 8.000 3rd Qu.:0.4188
## Max. :21.00 Max. :10.000 Max. :0.5250
## Number plotname
## Min. : 1.000 Length:104
## 1st Qu.: 1.000 Class :character
## Median : 2.000 Mode :character
## Mean : 3.298
## 3rd Qu.: 4.000
## Max. :16.000create species abundance matrix
alldata <- alldata %>% mutate(ForestType=as.factor(ForestType),Species=as.factor(Species),Elevation=as.factor(Elevation),plotname=as.factor(paste0(ForestType,No.Plot,Elevation)))
summary(alldata)## ForestType No.Plot Elevation Slope DistRiver
## Primary :55 Min. :1.00 500:32 Min. :0.500 Min. : 1.00
## Secondary:49 1st Qu.:1.00 600:28 1st Qu.:1.600 1st Qu.: 6.80
## Median :2.00 800:22 Median :4.000 Median :11.00
## Mean :2.01 900:22 Mean :3.945 Mean :10.35
## 3rd Qu.:3.00 3rd Qu.:6.000 3rd Qu.:13.00
## Max. :3.00 Max. :8.200 Max. :21.00
##
## CoordinateX CoordinateY Species
## Min. : 1.000 Min. :0.0500 Rhododendron_decorum :19
## 1st Qu.: 2.250 1st Qu.:0.1000 Itoa_orientalis :15
## Median : 4.000 Median :0.2125 Rhododendron_leptothrium:13
## Mean : 4.894 Mean :0.2719 Ficus_heteromorpha :10
## 3rd Qu.: 8.000 3rd Qu.:0.4188 Illicium_macranthum : 6
## Max. :10.000 Max. :0.5250 Manglietia_insignis : 5
## (Other) :36
## Number plotname
## Min. : 1.000 Secondary1500: 8
## 1st Qu.: 1.000 Primary2600 : 6
## Median : 2.000 Primary3600 : 6
## Mean : 3.298 Primary1500 : 5
## 3rd Qu.: 4.000 Primary2900 : 5
## Max. :16.000 Primary3500 : 5
## (Other) :69length(unique(alldata$Species))## [1] 22spabun<-as.data.frame(matrix(0,nrow=length(unique(alldata$plotname)),ncol=length(unique(alldata$Species))))
rownames(spabun)<-unique(alldata$plotname)
colnames(spabun)<-unique(alldata$Species)
summary(spabun)## Rhododendron_leptothrium Rhododendron_decorum Ficus_heteromorpha
## Min. :0 Min. :0 Min. :0
## 1st Qu.:0 1st Qu.:0 1st Qu.:0
## Median :0 Median :0 Median :0
## Mean :0 Mean :0 Mean :0
## 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0
## Max. :0 Max. :0 Max. :0
## Osbeckia_mepalensis Vaccinium_duclouxii Itoa_orientalis Illicium_macranthum
## Min. :0 Min. :0 Min. :0 Min. :0
## 1st Qu.:0 1st Qu.:0 1st Qu.:0 1st Qu.:0
## Median :0 Median :0 Median :0 Median :0
## Mean :0 Mean :0 Mean :0 Mean :0
## 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0
## Max. :0 Max. :0 Max. :0 Max. :0
## Manglietia_insignis Beilschmiedia_robusta Sterculia_nobilis
## Min. :0 Min. :0 Min. :0
## 1st Qu.:0 1st Qu.:0 1st Qu.:0
## Median :0 Median :0 Median :0
## Mean :0 Mean :0 Mean :0
## 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0
## Max. :0 Max. :0 Max. :0
## Saurauia_napaulensis Padus_napaulensis Craibiodendron_yunnanense
## Min. :0 Min. :0 Min. :0
## 1st Qu.:0 1st Qu.:0 1st Qu.:0
## Median :0 Median :0 Median :0
## Mean :0 Mean :0 Mean :0
## 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0
## Max. :0 Max. :0 Max. :0
## Rhododendron_irroratum Ilex_corallina Lithocarpus_hancei Schefflera_fengii
## Min. :0 Min. :0 Min. :0 Min. :0
## 1st Qu.:0 1st Qu.:0 1st Qu.:0 1st Qu.:0
## Median :0 Median :0 Median :0 Median :0
## Mean :0 Mean :0 Mean :0 Mean :0
## 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0
## Max. :0 Max. :0 Max. :0 Max. :0
## Skimmia_arborescens Litsea_chunii Claoxylon_khasianum Schefflera_shweliensis
## Min. :0 Min. :0 Min. :0 Min. :0
## 1st Qu.:0 1st Qu.:0 1st Qu.:0 1st Qu.:0
## Median :0 Median :0 Median :0 Median :0
## Mean :0 Mean :0 Mean :0 Mean :0
## 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0
## Max. :0 Max. :0 Max. :0 Max. :0
## Camellia_forrestii
## Min. :0
## 1st Qu.:0
## Median :0
## Mean :0
## 3rd Qu.:0
## Max. :0for(i in 1:nrow(alldata)){
spabun[as.character(alldata$plotname[i]),as.character(alldata$Species[i])]<-alldata$Number[i]
}
# write.csv(spabun,"SpeciesAbundance.csv")We altogether have 24 sites, and 22 species. Thus the species abundance matrix consists 24 rows and 22 columns.
look at trait
library(tidyverse)
trait <- read_csv("dummy_trait.csv")
summary(trait)## sp LMA LL Amass
## Length:77 Min. : 17.38 Min. : 5.38 Min. : 11.35
## Class :character 1st Qu.: 58.88 1st Qu.: 21.88 1st Qu.: 54.83
## Mode :character Median : 97.72 Median : 40.30 Median : 97.93
## Mean :132.88 Mean : 59.03 Mean :140.15
## 3rd Qu.:158.49 3rd Qu.: 81.85 3rd Qu.:174.55
## Max. :602.56 Max. :267.95 Max. :985.05
## Rmass Nmass Pmass WD
## Min. : 1.29 Min. : 0.300 Min. :0.0100 Min. :0.3100
## 1st Qu.: 7.23 1st Qu.: 0.960 1st Qu.:0.0500 1st Qu.:0.4400
## Median :13.94 Median : 1.530 Median :0.1100 Median :0.4900
## Mean :18.02 Mean : 2.744 Mean :0.1481 Mean :0.5018
## 3rd Qu.:23.80 3rd Qu.: 4.220 3rd Qu.:0.1600 3rd Qu.:0.5400
## Max. :67.58 Max. :22.690 Max. :0.7400 Max. :0.8600
## SM
## Min. :0.060
## 1st Qu.:0.550
## Median :1.350
## Mean :1.933
## 3rd Qu.:2.790
## Max. :7.980Trait columns meaning: LMA stands for Leaf mass per area (g m^2), LL stands for Leaf lifespans (longevity, months), Amass stands for Maximum photosynthetic rates per unit mass (nnmol g^-1 s^-1), Rmass stands for Dark respiration rates per unit mass (nnmol g^-1 s^-1), Nmass stands for Leaf nitrogen per unit mass (%), Pmass stands for Leaf phosphorus per unit mass (%), WD stands for Wood density (g cm^-3), SM stands for Seed dry mass (mg).
trait %>%
pivot_longer(LMA:SM, names_to = "trait") %>%
ggplot() + aes(x=value) +
geom_histogram(position = "identity") +
facet_wrap(~ trait, scale = "free")
# transforamtion
trait2 <- trait |>
mutate(logLMA = log(LMA),
logLL = log(LL),
logAmass = log(Amass),
logRmass = log(Rmass),
logNmass = log(Nmass),
logPmass = log(Pmass),
logSM = log(SM)) |>
dplyr::select(sp, logLMA, logLL, logAmass, logRmass, logNmass, logPmass, WD, logSM)
trait2 %>%
pivot_longer(logLMA:logSM, names_to = "trait") %>%
ggplot() + aes(x=value) +
geom_histogram(position = "identity") +
facet_wrap(~ trait, scale = "free")
# much normal nowDiversity indices calculating
diversity <- data.frame(plot=rownames(spabun))
diversity$shannon<-vegan::diversity(spabun, index = "shannon")
diversity$No.ind <- apply(spabun,1,sum)library(picante)
library(FD)
phylo <- read.tree("dummy_tree.newick")
plot(phylo)
diversity[,c("PD","richness")] <- pd(spabun, phylo)
diversity[,"MNTD"] <-mntd(spabun, cophenetic(phylo))
diversity[,"MPD"] <- mpd(spabun, cophenetic(phylo))trait_mat0 <- as.matrix(trait2[, -1])
rownames(trait_mat0) <- trait2$sp
# scale
trait_mat <- scale(trait_mat0)
tmp <- trait2 %>%
filter(sp %in% colnames(spabun))
res_fd <- dbFD(trait_mat[colnames(spabun), ], spabun)## FRic: Dimensionality reduction was required. The last 6 PCoA axes (out of 8 in total) were removed.
## FRic: Quality of the reduced-space representation = 0.8160306diversity$FRic <- res_fd$FRic
diversity$FDiv <- res_fd$FDiv
diversity$FDis <- res_fd$FDis
diversity$RaoQ <- res_fd$RaoQ
diversity$FEVe <- res_fd$FEveDiversity pattern across elevation and forest types
Now I have calculated all
library(stringr)
diversity$Elevation<-as.integer(str_sub(diversity$plot,start = str_length(diversity$plot)-2, end=-1))
diversity$forest <- as.factor(str_sub(diversity$plot,start = 1, end=str_length(diversity$plot)-4))
divcor <- cor(diversity[,2:12], method="pearson")
corrplot::corrplot(divcor)
Didn’t remove anything
diversity_long <- diversity %>% pivot_longer(shannon:RaoQ, names_to = "diversityindices") %>% mutate(diversityindices=factor(diversityindices, levels=c("No.ind","richness","shannon","PD","MPD","MNTD","FRic","FEVe","FDiv","FDis","RaoQ")))
ggplot(diversity_long, aes(x = Elevation, y = value, colour = forest)) + geom_point() + geom_smooth(method = lm) +facet_wrap(~diversityindices,scales = "free",ncol=3) + theme_classic()
library(lme4)
library(lmerTest)
library(performance)
ggplot(diversity_long,aes(x=value))+geom_histogram()+facet_wrap(~diversityindices,scale="free",ncol=3)
diversity %>% mutate(MNTD=log(MNTD),FRic=log(FRic),FEVe=(FEVe)^2,FDiv=(FDiv)^2,FDis=log(FDis),RaoQ=log(RaoQ)) %>%pivot_longer(shannon:RaoQ, names_to = "diversityindices") %>% mutate(diversityindices=factor(diversityindices, levels=c("No.ind","richness","shannon","PD","MPD","MNTD","FRic","FEVe","FDiv","FDis","RaoQ"))) %>% ggplot()+ aes(x = value) + geom_histogram() +facet_wrap(~diversityindices,scales = "free",ncol=3) + theme_classic()
#good enough now
diversity.trans<-diversity %>% mutate(MNTD=log(MNTD),FRic=log(FRic),FEVe=(FEVe)^2,FDiv=(FDiv)^2,FDis=log(FDis),RaoQ=log(RaoQ))
diversity.trans %>% pivot_longer(shannon:RaoQ, names_to = "diversityindices") %>% mutate(diversityindices=factor(diversityindices, levels=c("No.ind","richness","shannon","PD","MPD","MNTD","FRic","FEVe","FDiv","FDis","RaoQ"))) %>%ggplot()+aes(x = Elevation, y = value, colour = forest) + geom_point() + geom_smooth(method = lm) +facet_wrap(~diversityindices,scales = "free",ncol=3) + theme_classic()
diversity.trans$Elevation.sq <- (diversity.trans$Elevation)^2
##No. ind
indlm1 <- lm(No.ind~Elevation*forest, data=diversity.trans)
summary(indlm1) #signidicant##
## Call:
## lm(formula = No.ind ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.083 -1.371 -0.725 1.933 5.133
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -13.916667 3.789259 -3.673 0.00151 **
## Elevation 0.045000 0.005280 8.522 4.33e-08 ***
## forestSecondary 10.450000 5.358822 1.950 0.06533 .
## Elevation:forestSecondary -0.024333 0.007467 -3.259 0.00393 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.892 on 20 degrees of freedom
## Multiple R-squared: 0.8562, Adjusted R-squared: 0.8346
## F-statistic: 39.68 on 3 and 20 DF, p-value: 1.308e-08check_model(indlm1) # wrong!!!
indlm2 <- glm(No.ind~Elevation*forest, data=diversity.trans, family = "poisson")
summary(indlm2)##
## Call:
## glm(formula = No.ind ~ Elevation * forest, family = "poisson",
## data = diversity.trans)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.3143 -0.4419 -0.2295 0.3604 1.6520
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.9007468 0.3638483 2.476 0.0133 *
## Elevation 0.0026833 0.0004677 5.738 9.59e-09 ***
## forestSecondary 0.1031159 0.5675999 0.182 0.8558
## Elevation:forestSecondary -0.0007573 0.0007384 -1.026 0.3051
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 77.532 on 23 degrees of freedom
## Residual deviance: 12.014 on 20 degrees of freedom
## AIC: 125.53
##
## Number of Fisher Scoring iterations: 4check_model(indlm2)
indlm3 <- glm(No.ind~Elevation, data=diversity.trans, family = "poisson")
summary(indlm3)##
## Call:
## glm(formula = No.ind ~ Elevation, family = "poisson", data = diversity.trans)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.9752 -0.8117 -0.1612 0.9443 2.2533
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.9198876 0.2790661 3.296 0.00098 ***
## Elevation 0.0023857 0.0003615 6.600 4.11e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 77.532 on 23 degrees of freedom
## Residual deviance: 31.422 on 22 degrees of freedom
## AIC: 140.94
##
## Number of Fisher Scoring iterations: 4check_model(indlm3)
indlm4 <- update(indlm2,.~.+Elevation.sq*forest)
summary(indlm4) #not siignificant##
## Call:
## glm(formula = No.ind ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, family = "poisson", data = diversity.trans)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.2850 -0.4980 -0.1844 0.3629 1.5247
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.266e+00 2.316e+00 0.547 0.585
## Elevation 1.605e-03 6.766e-03 0.237 0.813
## forestSecondary 8.432e-01 3.670e+00 0.230 0.818
## Elevation.sq 7.571e-07 4.738e-06 0.160 0.873
## Elevation:forestSecondary -2.963e-03 1.077e-02 -0.275 0.783
## forestSecondary:Elevation.sq 1.560e-06 7.569e-06 0.206 0.837
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 77.532 on 23 degrees of freedom
## Residual deviance: 11.835 on 18 degrees of freedom
## AIC: 129.35
##
## Number of Fisher Scoring iterations: 4AIC(indlm2,indlm3) # Although AIC lower, indlm3 diagnostic better## df AIC
## indlm2 4 125.5309
## indlm3 2 140.9388##richness
richnesslm1 <- lm(richness~Elevation*forest, data=diversity.trans)
summary(richnesslm1)##
## Call:
## lm(formula = richness ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.7167 -0.4500 -0.2333 0.3750 2.6500
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.216667 1.315664 4.725 0.00013 ***
## Elevation -0.002333 0.001833 -1.273 0.21771
## forestSecondary 2.300000 1.860630 1.236 0.23073
## Elevation:forestSecondary -0.004000 0.002593 -1.543 0.13856
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.004 on 20 degrees of freedom
## Multiple R-squared: 0.4292, Adjusted R-squared: 0.3436
## F-statistic: 5.014 on 3 and 20 DF, p-value: 0.009405#not significant
check_model(richnesslm1)# wrong!!!
richnesslm2 <- glm(richness~Elevation*forest, data=diversity.trans, family = "poisson")
summary(richnesslm2)##
## Call:
## glm(formula = richness ~ Elevation * forest, family = "poisson",
## data = diversity.trans)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -0.8126 -0.2116 -0.1141 0.2024 1.0306
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.8761749 0.6029419 3.112 0.00186 **
## Elevation -0.0005100 0.0008551 -0.596 0.55089
## forestSecondary 0.6039792 0.8713886 0.693 0.48823
## Elevation:forestSecondary -0.0010673 0.0012608 -0.846 0.39729
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 7.7100 on 23 degrees of freedom
## Residual deviance: 4.0354 on 20 degrees of freedom
## AIC: 91.456
##
## Number of Fisher Scoring iterations: 4#not significant
richnesslm3 <-update(richnesslm2,.~.+Elevation.sq*forest)
summary(richnesslm3)#not significant##
## Call:
## glm(formula = richness ~ Elevation + forest + Elevation.sq +
## Elevation:forest + forest:Elevation.sq, family = "poisson",
## data = diversity.trans)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -0.52760 -0.30006 0.00513 0.23344 0.77838
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.298e-01 4.213e+00 0.031 0.975
## Elevation 4.763e-03 1.261e-02 0.378 0.706
## forestSecondary 6.235e+00 6.133e+00 1.017 0.309
## Elevation.sq -3.779e-06 9.014e-06 -0.419 0.675
## Elevation:forestSecondary -1.818e-02 1.847e-02 -0.984 0.325
## forestSecondary:Elevation.sq 1.231e-05 1.325e-05 0.930 0.353
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 7.7100 on 23 degrees of freedom
## Residual deviance: 3.0798 on 18 degrees of freedom
## AIC: 94.5
##
## Number of Fisher Scoring iterations: 4##shannon
shannonlm1 <- lm(shannon~Elevation*forest,data=diversity.trans)
summary(shannonlm1)#not significant##
## Call:
## lm(formula = shannon ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.45016 -0.15297 -0.01576 0.19303 0.43486
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.7561507 0.3216850 5.459 2.41e-05 ***
## Elevation -0.0006424 0.0004483 -1.433 0.1673
## forestSecondary 0.7578343 0.4549313 1.666 0.1113
## Elevation:forestSecondary -0.0012969 0.0006339 -2.046 0.0542 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2455 on 20 degrees of freedom
## Multiple R-squared: 0.535, Adjusted R-squared: 0.4652
## F-statistic: 7.67 on 3 and 20 DF, p-value: 0.001328shannonlm2 <-update(shannonlm1,.~.+Elevation.sq*forest)
summary(shannonlm2)#significant##
## Call:
## lm(formula = shannon ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.2642 -0.1451 -0.0284 0.1365 0.3963
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4.046e-01 1.824e+00 -0.222 0.82697
## Elevation 5.863e-03 5.447e-03 1.076 0.29594
## forestSecondary 8.684e+00 2.580e+00 3.366 0.00344 **
## Elevation.sq -4.647e-06 3.882e-06 -1.197 0.24678
## Elevation:forestSecondary -2.516e-02 7.703e-03 -3.266 0.00429 **
## forestSecondary:Elevation.sq 1.705e-05 5.489e-06 3.105 0.00611 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2017 on 18 degrees of freedom
## Multiple R-squared: 0.7176, Adjusted R-squared: 0.6391
## F-statistic: 9.147 on 5 and 18 DF, p-value: 0.0001817check_model(shannonlm2) #fine
##PD
PDlm1 <- lm(PD~Elevation*forest, data=diversity.trans)
summary(PDlm1)##
## Call:
## lm(formula = PD ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.03513 -0.27121 -0.08996 0.32566 1.33305
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.5841856 0.8017674 5.718 1.35e-05 ***
## Elevation -0.0029754 0.0011172 -2.663 0.0149 *
## forestSecondary 0.5429924 1.1338703 0.479 0.6372
## Elevation:forestSecondary -0.0007405 0.0015800 -0.469 0.6444
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.6119 on 20 degrees of freedom
## Multiple R-squared: 0.476, Adjusted R-squared: 0.3973
## F-statistic: 6.055 on 3 and 20 DF, p-value: 0.004179PDlm2 <- lm(PD~Elevation, data=diversity.trans)
summary(PDlm2)##
## Call:
## lm(formula = PD ~ Elevation, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.05985 -0.30268 -0.09394 0.36579 1.41941
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.8556818 0.5436424 8.932 9.04e-09 ***
## Elevation -0.0033456 0.0007575 -4.416 0.000218 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5868 on 22 degrees of freedom
## Multiple R-squared: 0.4699, Adjusted R-squared: 0.4458
## F-statistic: 19.5 on 1 and 22 DF, p-value: 0.0002182PDlm3 <- update(PDlm1,.~.+Elevation.sq*forest)
summary(PDlm3) #not significant##
## Call:
## lm(formula = PD ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.0358 -0.3317 -0.0536 0.3694 1.2540
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 7.182e+00 5.772e+00 1.244 0.229
## Elevation -1.080e-02 1.723e-02 -0.627 0.539
## forestSecondary 3.962e-01 8.162e+00 0.049 0.962
## Elevation.sq 5.587e-06 1.228e-05 0.455 0.655
## Elevation:forestSecondary -2.986e-04 2.437e-02 -0.012 0.990
## forestSecondary:Elevation.sq -3.157e-07 1.737e-05 -0.018 0.986
##
## Residual standard error: 0.6381 on 18 degrees of freedom
## Multiple R-squared: 0.4871, Adjusted R-squared: 0.3446
## F-statistic: 3.419 on 5 and 18 DF, p-value: 0.024AIC(PDlm1,PDlm2)## df AIC
## PDlm1 5 50.15914
## PDlm2 3 46.43283anova(PDlm1,PDlm2)## Analysis of Variance Table
##
## Model 1: PD ~ Elevation * forest
## Model 2: PD ~ Elevation
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 20 7.4893
## 2 22 7.5752 -2 -0.085895 0.1147 0.8922#just keep elevation
##MPD
MPDlm1 <- lm(MPD~Elevation*forest,data=diversity.trans)
summary(MPDlm1)#not significant##
## Call:
## lm(formula = MPD ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.95745 -0.14704 0.02182 0.11525 0.54686
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.9718876 0.4656388 4.235 0.000406 ***
## Elevation -0.0012258 0.0006489 -1.889 0.073464 .
## forestSecondary -0.1566694 0.6585127 -0.238 0.814368
## Elevation:forestSecondary 0.0003193 0.0009176 0.348 0.731529
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3554 on 20 degrees of freedom
## Multiple R-squared: 0.2228, Adjusted R-squared: 0.1062
## F-statistic: 1.911 on 3 and 20 DF, p-value: 0.1604MPDlm2 <- update(MPDlm1,.~.+Elevation.sq*forest)
summary(MPDlm2)#not significant##
## Call:
## lm(formula = MPD ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.03991 -0.10106 0.01737 0.12620 0.52058
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.601e+00 3.274e+00 1.406 0.177
## Elevation -9.142e-03 9.774e-03 -0.935 0.362
## forestSecondary -5.342e+00 4.630e+00 -1.154 0.264
## Elevation.sq 5.654e-06 6.966e-06 0.812 0.428
## Elevation:forestSecondary 1.593e-02 1.382e-02 1.153 0.264
## forestSecondary:Elevation.sq -1.115e-05 9.851e-06 -1.132 0.272
##
## Residual standard error: 0.3619 on 18 degrees of freedom
## Multiple R-squared: 0.2744, Adjusted R-squared: 0.0729
## F-statistic: 1.362 on 5 and 18 DF, p-value: 0.2846##MNTD
MNTDlm1 <- lm(MNTD~Elevation*forest,data=diversity.trans)
summary(MNTDlm1)#not significant##
## Call:
## lm(formula = MNTD ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.70494 -0.15764 0.04962 0.19746 0.91350
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.6952460 0.7429117 0.936 0.3605
## Elevation -0.0018248 0.0010352 -1.763 0.0932 .
## forestSecondary -0.2976132 1.0506358 -0.283 0.7799
## Elevation:forestSecondary 0.0006542 0.0014640 0.447 0.6598
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.567 on 20 degrees of freedom
## Multiple R-squared: 0.1957, Adjusted R-squared: 0.07502
## F-statistic: 1.622 on 3 and 20 DF, p-value: 0.2159MNTDlm2 <-update(MNTDlm1,.~.+Elevation.sq*forest)
summary(MNTDlm2)#not significant##
## Call:
## lm(formula = MNTD ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.81234 -0.12855 0.02618 0.20882 0.80610
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.823e+00 5.206e+00 1.118 0.278
## Elevation -1.726e-02 1.554e-02 -1.111 0.281
## forestSecondary -8.755e+00 7.362e+00 -1.189 0.250
## Elevation.sq 1.103e-05 1.108e-05 0.995 0.333
## Elevation:forestSecondary 2.612e-02 2.198e-02 1.188 0.250
## forestSecondary:Elevation.sq -1.819e-05 1.567e-05 -1.161 0.261
##
## Residual standard error: 0.5756 on 18 degrees of freedom
## Multiple R-squared: 0.254, Adjusted R-squared: 0.04684
## F-statistic: 1.226 on 5 and 18 DF, p-value: 0.3374##FRic
FRiclm1 <- lm(FRic~Elevation*forest,data=diversity.trans)
summary(FRiclm1) #not significant##
## Call:
## lm(formula = FRic ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.8991 -0.1224 0.1448 0.5541 1.0231
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.5044259 1.1066559 1.359 0.189
## Elevation -0.0002271 0.0015421 -0.147 0.884
## forestSecondary 0.9192560 1.5650477 0.587 0.564
## Elevation:forestSecondary -0.0011251 0.0021808 -0.516 0.612
##
## Residual standard error: 0.8446 on 20 degrees of freedom
## Multiple R-squared: 0.04473, Adjusted R-squared: -0.09856
## F-statistic: 0.3121 on 3 and 20 DF, p-value: 0.8164FRiclm2 <- update(FRiclm1, .~.+Elevation.sq*forest)
summary(FRiclm2) # significant##
## Call:
## lm(formula = FRic ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.4962 -0.3510 0.1637 0.4185 1.0458
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.099e+01 7.031e+00 -1.563 0.1355
## Elevation 3.738e-02 2.099e-02 1.781 0.0918 .
## forestSecondary 2.416e+01 9.943e+00 2.430 0.0258 *
## Elevation.sq -2.686e-05 1.496e-05 -1.796 0.0893 .
## Elevation:forestSecondary -7.110e-02 2.969e-02 -2.395 0.0277 *
## forestSecondary:Elevation.sq 4.998e-05 2.116e-05 2.362 0.0296 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7774 on 18 degrees of freedom
## Multiple R-squared: 0.2718, Adjusted R-squared: 0.06949
## F-statistic: 1.344 on 5 and 18 DF, p-value: 0.2912check_model(FRiclm2) # diagnostic not good
##FEVe
FEVelm1 <- lm(FEVe~Elevation*forest,data=diversity.trans)
summary(FEVelm1)##
## Call:
## lm(formula = FEVe ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.29327 -0.09582 0.01052 0.12074 0.29076
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.7660142 0.2237837 3.423 0.00269 **
## Elevation -0.0003705 0.0003118 -1.188 0.24868
## forestSecondary 0.2298707 0.3164779 0.726 0.47605
## Elevation:forestSecondary -0.0004105 0.0004410 -0.931 0.36299
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1708 on 20 degrees of freedom
## Multiple R-squared: 0.2949, Adjusted R-squared: 0.1892
## F-statistic: 2.789 on 3 and 20 DF, p-value: 0.06712#not significant
FEVelm2<-update(FEVelm1, .~.+Elevation.sq*forest)
summary(FEVelm2)# not significant##
## Call:
## lm(formula = FEVe ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.253418 -0.119571 0.005495 0.110520 0.283927
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4.693e-01 1.601e+00 -0.293 0.773
## Elevation 3.349e-03 4.779e-03 0.701 0.492
## forestSecondary 1.254e+00 2.264e+00 0.554 0.587
## Elevation.sq -2.657e-06 3.406e-06 -0.780 0.446
## Elevation:forestSecondary -3.493e-03 6.759e-03 -0.517 0.612
## forestSecondary:Elevation.sq 2.201e-06 4.817e-06 0.457 0.653
##
## Residual standard error: 0.177 on 18 degrees of freedom
## Multiple R-squared: 0.3186, Adjusted R-squared: 0.1294
## F-statistic: 1.683 on 5 and 18 DF, p-value: 0.1896##FDiv
FDivlm1 <- lm(FDiv~Elevation*forest,data=diversity.trans)
summary(FDivlm1)# significant##
## Call:
## lm(formula = FDiv ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.30058 -0.08177 0.04738 0.10321 0.15507
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.4506895 0.1770097 2.546 0.0192 *
## Elevation 0.0002208 0.0002467 0.895 0.3814
## forestSecondary 0.5176315 0.2503295 2.068 0.0518 .
## Elevation:forestSecondary -0.0007461 0.0003488 -2.139 0.0450 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1351 on 20 degrees of freedom
## Multiple R-squared: 0.2109, Adjusted R-squared: 0.09249
## F-statistic: 1.781 on 3 and 20 DF, p-value: 0.1831check_model(FDivlm1) # not good enough but make sense
FDivlm2 <- update(FDivlm1,.~.+Elevation.sq*forest)
summary(FDivlm2) # not significant##
## Call:
## lm(formula = FDiv ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.25571 -0.07467 0.01305 0.09138 0.18628
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -9.403e-01 1.235e+00 -0.761 0.456
## Elevation 4.409e-03 3.688e-03 1.195 0.247
## forestSecondary 1.265e+00 1.747e+00 0.724 0.478
## Elevation.sq -2.991e-06 2.628e-06 -1.138 0.270
## Elevation:forestSecondary -2.996e-03 5.216e-03 -0.574 0.573
## forestSecondary:Elevation.sq 1.607e-06 3.717e-06 0.432 0.671
##
## Residual standard error: 0.1366 on 18 degrees of freedom
## Multiple R-squared: 0.2743, Adjusted R-squared: 0.07268
## F-statistic: 1.361 on 5 and 18 DF, p-value: 0.285##FDis
FDislm1 <- lm(FDis~Elevation*forest,data=diversity.trans)
summary(FDislm1)# significanyt##
## Call:
## lm(formula = FDis ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.57008 -0.17965 0.04284 0.18087 0.30536
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.7423284 0.3379152 2.197 0.0400 *
## Elevation -0.0003485 0.0004709 -0.740 0.4679
## forestSecondary 1.0471963 0.4778842 2.191 0.0404 *
## Elevation:forestSecondary -0.0014660 0.0006659 -2.202 0.0396 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2579 on 20 degrees of freedom
## Multiple R-squared: 0.4356, Adjusted R-squared: 0.351
## F-statistic: 5.146 on 3 and 20 DF, p-value: 0.008458FDislm2 <-update(FDislm1, .~.+Elevation.sq*forest)
summary(FDislm2) # very very significant##
## Call:
## lm(formula = FDis ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.39544 -0.09303 0.03150 0.09719 0.32101
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4.229e+00 1.726e+00 -2.450 0.024762 *
## Elevation 1.462e-02 5.154e-03 2.836 0.010954 *
## forestSecondary 1.143e+01 2.441e+00 4.683 0.000185 ***
## Elevation.sq -1.069e-05 3.673e-06 -2.910 0.009334 **
## Elevation:forestSecondary -3.273e-02 7.289e-03 -4.490 0.000283 ***
## forestSecondary:Elevation.sq 2.233e-05 5.195e-06 4.299 0.000432 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1909 on 18 degrees of freedom
## Multiple R-squared: 0.7218, Adjusted R-squared: 0.6445
## F-statistic: 9.34 on 5 and 18 DF, p-value: 0.0001599check_model(FDislm2) #good
##RaoQ
RaoQlm1 <-lm(RaoQ~Elevation*forest,data=diversity.trans)
summary(RaoQlm1) # not significance##
## Call:
## lm(formula = RaoQ ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.1537 -0.2838 0.1330 0.3700 0.5420
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.8620844 0.6464953 2.880 0.00925 **
## Elevation -0.0009568 0.0009009 -1.062 0.30086
## forestSecondary 1.2459579 0.9142824 1.363 0.18810
## Elevation:forestSecondary -0.0015943 0.0012740 -1.251 0.22523
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4934 on 20 degrees of freedom
## Multiple R-squared: 0.3235, Adjusted R-squared: 0.222
## F-statistic: 3.188 on 3 and 20 DF, p-value: 0.04599RaoQlm2 <- update(RaoQlm1,.~.+Elevation.sq*forest)
summary(RaoQlm2) # significant##
## Call:
## lm(formula = RaoQ ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.85525 -0.20003 0.08947 0.22197 0.52818
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6.822e+00 3.605e+00 -1.892 0.07462 .
## Elevation 2.519e-02 1.076e-02 2.340 0.03099 *
## forestSecondary 1.918e+01 5.098e+00 3.763 0.00142 **
## Elevation.sq -1.867e-05 7.670e-06 -2.435 0.02553 *
## Elevation:forestSecondary -5.559e-02 1.522e-02 -3.652 0.00182 **
## forestSecondary:Elevation.sq 3.857e-05 1.085e-05 3.556 0.00226 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3985 on 18 degrees of freedom
## Multiple R-squared: 0.6028, Adjusted R-squared: 0.4925
## F-statistic: 5.463 on 5 and 18 DF, p-value: 0.003131check_model(RaoQlm2)# not good
Significant models: indlm3, shannonlm2, PDlm2, FDivlm1, FDislm2, last two transformed
indlm3$call## glm(formula = No.ind ~ Elevation, family = "poisson", data = diversity.trans)check_overdispersion(indlm3) #no overdispersion## # Overdispersion test
##
## dispersion ratio = 1.426
## Pearson's Chi-Squared = 31.380
## p-value = 0.089indnew <- data.frame(Elevation=rep(seq(500,900,10),2),forest=rep(c("Primary","Secondary"),each=41))
indnew$Elevation.sq <-(indnew$Elevation)^2
preds <- predict(indlm3, newdata=indnew,type= 'link', se.fit=T)
indnew$No.ind <- exp(preds$fit)
indnew$indup <- exp(preds$fit+preds$se.fit*1.96)
indnew$indlw <- exp(preds$fit-preds$se.fit*1.96)
indp <- ggplot(diversity,aes(x=Elevation,y=No.ind,color=forest))+geom_point()+geom_smooth(data=indnew,aes(x=Elevation, y=No.ind, ymin=indlw, ymax=indup),stat='identity')+theme_classic()
shannonlm2$call## lm(formula = shannon ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)shannonnew <- data.frame(Elevation=rep(seq(500,900,10),2),forest=rep(c("Primary","Secondary"),each=41))
shannonnew$Elevation.sq <-(shannonnew$Elevation)^2
preds <- predict(shannonlm2, newdata=shannonnew, se.fit=T)
shannonnew$shannon <- preds$fit
shannonnew$shannonup <- preds$fit+preds$se.fit*1.96
shannonnew$shannonlw <- preds$fit-preds$se.fit*1.96
shannonp <- ggplot(diversity,aes(x=Elevation,y=shannon,color=forest))+geom_point()+geom_smooth(data=shannonnew,aes(x=Elevation, y=shannon, ymin=shannonlw, ymax=shannonup),stat='identity')+theme_classic()
PDlm2$call## lm(formula = PD ~ Elevation, data = diversity.trans)summary(PDlm2)##
## Call:
## lm(formula = PD ~ Elevation, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.05985 -0.30268 -0.09394 0.36579 1.41941
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.8556818 0.5436424 8.932 9.04e-09 ***
## Elevation -0.0033456 0.0007575 -4.416 0.000218 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5868 on 22 degrees of freedom
## Multiple R-squared: 0.4699, Adjusted R-squared: 0.4458
## F-statistic: 19.5 on 1 and 22 DF, p-value: 0.0002182PDnew <- data.frame(Elevation=rep(seq(500,900,10),2),forest=rep(c("Primary","Secondary"),each=41))
PDnew$Elevation.sq <-(PDnew$Elevation)^2
preds <- predict(PDlm2, newdata=PDnew, se.fit=T)
PDnew$PD <- preds$fit
PDnew$PDup <- preds$fit+preds$se.fit*1.96
PDnew$PDlw <- preds$fit-preds$se.fit*1.96
PDp <- ggplot(diversity,aes(x=Elevation,y=PD,color=forest))+geom_point()+geom_smooth(data=PDnew,aes(x=Elevation, y=PD, ymin=PDlw, ymax=PDup),stat='identity')+theme_classic()
FDivlm1$call## lm(formula = FDiv ~ Elevation * forest, data = diversity.trans)summary(FDivlm1)##
## Call:
## lm(formula = FDiv ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.30058 -0.08177 0.04738 0.10321 0.15507
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.4506895 0.1770097 2.546 0.0192 *
## Elevation 0.0002208 0.0002467 0.895 0.3814
## forestSecondary 0.5176315 0.2503295 2.068 0.0518 .
## Elevation:forestSecondary -0.0007461 0.0003488 -2.139 0.0450 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1351 on 20 degrees of freedom
## Multiple R-squared: 0.2109, Adjusted R-squared: 0.09249
## F-statistic: 1.781 on 3 and 20 DF, p-value: 0.1831FDivnew <- data.frame(Elevation=rep(seq(500,900,10),2),forest=rep(c("Primary","Secondary"),each=41))
FDivnew$Elevation.sq <-(FDivnew$Elevation)^2
preds <- predict(FDivlm1, newdata=FDivnew, se.fit=T)
FDivnew$FDiv <- sqrt(preds$fit)
FDivnew$FDivup <- sqrt(preds$fit+preds$se.fit*1.96)
FDivnew$FDivlw <- sqrt(preds$fit-preds$se.fit*1.96)
FDivp <- ggplot(diversity,aes(x=Elevation,y=FDiv,color=forest))+geom_point()+geom_smooth(data=FDivnew,aes(x=Elevation, y=FDiv, ymin=FDivlw, ymax=FDivup),stat='identity')+theme_classic()
FDislm2$call## lm(formula = FDis ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)summary(FDislm2)##
## Call:
## lm(formula = FDis ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.39544 -0.09303 0.03150 0.09719 0.32101
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4.229e+00 1.726e+00 -2.450 0.024762 *
## Elevation 1.462e-02 5.154e-03 2.836 0.010954 *
## forestSecondary 1.143e+01 2.441e+00 4.683 0.000185 ***
## Elevation.sq -1.069e-05 3.673e-06 -2.910 0.009334 **
## Elevation:forestSecondary -3.273e-02 7.289e-03 -4.490 0.000283 ***
## forestSecondary:Elevation.sq 2.233e-05 5.195e-06 4.299 0.000432 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1909 on 18 degrees of freedom
## Multiple R-squared: 0.7218, Adjusted R-squared: 0.6445
## F-statistic: 9.34 on 5 and 18 DF, p-value: 0.0001599FDisnew <- data.frame(Elevation=rep(seq(500,900,10),2),forest=rep(c("Primary","Secondary"),each=41))
FDisnew$Elevation.sq <-(FDisnew$Elevation)^2
preds <- predict(FDislm2, newdata=FDisnew, se.fit=T)
FDisnew$FDis <- exp(preds$fit)
FDisnew$FDisup <- exp(preds$fit+preds$se.fit*1.96)
FDisnew$FDislw <- exp(preds$fit-preds$se.fit*1.96)
FDisp <- ggplot(diversity,aes(x=Elevation,y=FDis,color=forest))+geom_point()+geom_smooth(data=FDisnew,aes(x=Elevation, y=FDis, ymin=FDislw, ymax=FDisup),stat='identity')+theme_classic()
library(patchwork)
indp+shannonp+PDp+FDivp+FDisp+plot_layout(ncol = 2,guides="collect")
Functional composition change
library(vegan)
traitsite <- res_fd$CWM
head(traitsite)## logLMA logLL logAmass logRmass logNmass logPmass
## Primary2900 1.985436 1.6514668 -1.5221013 -1.485665 -1.611124 -1.5231640
## Primary3900 1.940231 1.6073298 -1.6568297 -1.455187 -1.456864 -1.5736016
## Secondary1900 1.446783 1.1666166 -0.8020430 -1.525994 -1.405698 -0.9752096
## Secondary2900 1.150931 0.9119735 -0.6031357 -1.355206 -1.159585 -0.7083924
## Secondary3900 1.500520 1.2193407 -0.9559353 -1.628596 -1.330057 -1.2839945
## Primary1900 2.190338 1.8213643 -1.8174421 -1.580867 -1.678590 -1.7362137
## WD logSM
## Primary2900 0.441978928 0.5686285
## Primary3900 0.001045691 0.6428736
## Secondary1900 1.064215148 0.7003100
## Secondary2900 0.842350343 0.7870687
## Secondary3900 0.833067538 0.5465483
## Primary1900 0.184690749 0.5906153traitsite$Elevation <- diversity$Elevation
traitsite$forest <- diversity$forest
# calculate trait distance
trait_mat0 <- as.matrix(trait2[, -1])
rownames(trait_mat0) <- trait2$sp
trait_mat0[1:5, 1:5]## logLMA logLL logAmass logRmass logNmass
## Acer_campbellii 3.684118 1.957274 6.892692 4.002047 1.8809906
## Actinodaphne_forrestii 4.236712 2.527327 5.006359 2.173615 0.4121097
## Alnus_nepalensis 4.743366 4.010419 4.341335 2.022871 0.5007753
## Anneslea_fragrans 4.190715 3.293241 5.162211 3.703522 1.4632554
## Beilschmiedia_robusta 3.614964 3.085573 5.722441 3.526655 1.7544037trait_mat <- scale(trait_mat0)
trait_dm <- as.matrix(dist(trait_mat))
pcaall <- prcomp(traitsite[, c("logLMA","logPmass","WD","logNmass")], scale = T)
summary(pcaall)## Importance of components:
## PC1 PC2 PC3 PC4
## Standard deviation 1.6894 0.9810 0.40875 0.12858
## Proportion of Variance 0.7135 0.2406 0.04177 0.00413
## Cumulative Proportion 0.7135 0.9541 0.99587 1.00000biplot(pcaall)
corrplot::corrplot(cor(trait2[,2:9]))
library(ggfortify)
autoplot(pcaall, data = traitsite, colour = 'Elevation',shape="forest",loadings = TRUE, loadings.colour = 'blue',loadings.label = TRUE, loadings.label.size = 3)
##RDA
env <- read_csv("env data(1).csv") %>% select(Site,Slope,DistRiver,soil_n,soil_p,soil_moist,Elevation) %>% mutate(Slope=1/Slope)
corrplot::corrplot(cor(env[,c(2:7)]))
trait.rda <- rda(traitsite[,1:8] ~ Slope+DistRiver+soil_n+soil_p+soil_moist, data=env,scale=T)
summary(trait.rda)##
## Call:
## rda(formula = traitsite[, 1:8] ~ Slope + DistRiver + soil_n + soil_p + soil_moist, data = env, scale = T)
##
## Partitioning of correlations:
## Inertia Proportion
## Total 8.000 1.0000
## Constrained 2.703 0.3379
## Unconstrained 5.297 0.6621
##
## Eigenvalues, and their contribution to the correlations
##
## Importance of components:
## RDA1 RDA2 RDA3 RDA4 RDA5 PC1 PC2
## Eigenvalue 1.9300 0.55963 0.1760 0.036197 0.0016794 3.7710 0.9066
## Proportion Explained 0.2412 0.06995 0.0220 0.004525 0.0002099 0.4714 0.1133
## Cumulative Proportion 0.2412 0.31120 0.3332 0.337724 0.3379339 0.8093 0.9226
## PC3 PC4 PC5 PC6 PC7 PC8
## Eigenvalue 0.32417 0.17539 0.0984 0.017610 0.0021691 0.0011705
## Proportion Explained 0.04052 0.02192 0.0123 0.002201 0.0002711 0.0001463
## Cumulative Proportion 0.96316 0.98508 0.9974 0.999583 0.9998537 1.0000000
##
## Accumulated constrained eigenvalues
## Importance of components:
## RDA1 RDA2 RDA3 RDA4 RDA5
## Eigenvalue 1.9300 0.5596 0.1760 0.03620 0.0016794
## Proportion Explained 0.7139 0.2070 0.0651 0.01339 0.0006212
## Cumulative Proportion 0.7139 0.9209 0.9860 0.99938 1.0000000
##
## Scaling 2 for species and site scores
## * Species are scaled proportional to eigenvalues
## * Sites are unscaled: weighted dispersion equal on all dimensions
## * General scaling constant of scores: 3.683023
##
##
## Species scores
##
## RDA1 RDA2 RDA3 RDA4 RDA5 PC1
## logLMA 0.76109 0.10836 -0.05665 0.040301 -0.010244 -1.0373
## logLL 0.75975 0.11457 -0.05377 0.056952 -0.005212 -1.0299
## logAmass -0.73113 -0.29358 -0.01416 -0.004176 0.014088 0.9416
## logRmass -0.69284 0.10450 -0.04093 0.081831 -0.043736 1.0021
## logNmass -0.68082 -0.01708 0.13628 -0.115469 -0.012084 1.0564
## logPmass -0.70270 -0.02616 -0.16650 0.159191 0.017454 0.9388
## WD 0.08344 -0.31776 -0.46680 -0.089825 -0.007834 -0.4699
## logSM 0.37025 -0.85146 0.16200 0.055222 -0.009600 -0.3851
##
##
## Site scores (weighted sums of species scores)
##
## RDA1 RDA2 RDA3 RDA4 RDA5 PC1
## Primary2900 1.233980 0.05467 -0.519813 -0.24327 -6.5039 -1.29999
## Primary3900 1.207463 0.25466 1.051976 0.47751 -7.5744 -0.04737
## Secondary1900 0.486935 -1.13089 -2.352917 -1.39101 3.5155 -0.50369
## Secondary2900 0.031632 -1.32015 -1.571822 -0.48336 4.5468 -0.48435
## Secondary3900 0.641295 -0.61755 -1.509279 -2.55247 3.2870 -0.29217
## Primary1900 1.571049 0.37571 0.359002 0.08622 -8.2777 -0.97494
## Secondary1600 0.221980 0.13356 1.619350 -3.60535 5.2876 -0.53411
## Secondary2600 -2.163666 -0.88396 3.393189 2.59004 -5.9402 0.50954
## Secondary3600 0.500784 2.12290 1.200711 -0.73403 10.7746 -0.43496
## Primary1600 -0.307541 -1.23573 0.001178 1.07846 3.2856 -0.39620
## Primary2600 -0.041258 -0.33569 1.965297 -1.31328 0.7281 0.61740
## Primary3600 -0.283479 -0.41958 0.931236 1.37347 -7.7930 0.18654
## Primary1500 -0.182050 1.68598 -2.971604 2.00103 -4.3857 -0.54514
## Primary2500 0.494289 1.65862 0.992416 3.09112 4.2738 -0.04527
## Primary3500 0.495089 2.09490 1.431710 -1.00241 15.0008 -0.04361
## Secondary1500 -1.183557 0.05789 0.251527 0.97077 -1.5737 1.30926
## Secondary2500 -1.975871 1.16052 -0.364507 -3.60748 -3.2899 0.59514
## Secondary3500 -4.307248 0.51056 1.105889 0.96371 -5.7733 2.45841
## Secondary1800 -0.191235 -1.59397 -1.969326 1.56695 8.8351 -0.25386
## Secondary2800 1.035947 -0.63964 -2.323364 -3.21173 2.1827 -0.32220
## Secondary3800 0.001795 -1.65000 -2.545703 0.98423 6.5876 0.45876
## Primary1800 0.790354 0.02749 1.551853 1.24679 -6.1950 0.38255
## Primary2800 1.337037 0.13537 0.380758 1.86419 -8.7035 -0.32017
## Primary3800 0.586274 -0.44564 -0.107759 -0.15012 -2.2948 -0.01956
##
##
## Site constraints (linear combinations of constraining variables)
##
## RDA1 RDA2 RDA3 RDA4 RDA5 PC1
## Primary2900 -0.60367 -0.59178 -0.53712 1.4930989 -1.63768 -1.29999
## Primary3900 1.10225 1.35943 0.17342 0.4550296 0.43850 -0.04737
## Secondary1900 -0.04255 0.05827 1.21642 0.1519446 0.69571 -0.50369
## Secondary2900 -0.67769 -1.40917 -2.14041 -0.5152921 1.52326 -0.48435
## Secondary3900 0.29124 0.36152 -0.72704 -0.0123036 -0.77557 -0.29217
## Primary1900 0.15357 0.14903 -0.12455 -1.4232908 -0.75958 -0.97494
## Secondary1600 -0.55043 -0.31086 0.96233 -1.6051563 -0.25621 -0.53411
## Secondary2600 -1.56998 -0.51907 0.89556 0.4548021 0.46655 0.50954
## Secondary3600 -0.16167 1.72883 -0.13149 0.5935009 0.65281 -0.43496
## Primary1600 -0.80481 -0.82919 1.42588 0.2669408 -0.49215 -0.39620
## Primary2600 0.77341 -0.51957 0.54220 -0.6723882 0.03780 0.61740
## Primary3600 -0.03167 -0.24091 0.86823 0.0592021 0.68517 0.18654
## Primary1500 -0.86892 1.35806 -0.70800 0.0129324 -0.02258 -0.54514
## Primary2500 0.35998 0.59774 -0.00142 -0.1297961 0.42925 -0.04527
## Primary3500 0.39877 0.17232 0.59514 0.5696806 -0.17479 -0.04361
## Secondary1500 0.62988 -0.62136 -0.17343 0.2409680 0.64127 1.30926
## Secondary2500 -1.12321 0.92455 -0.24831 -1.5151454 -0.03152 0.59514
## Secondary3500 -0.90471 0.37067 -0.49886 0.6326066 -0.75233 2.45841
## Secondary1800 -0.47117 -0.28803 -0.44118 0.9820175 -0.16854 -0.25386
## Secondary2800 0.70309 -0.48469 -0.19620 0.2202266 0.74500 -0.32220
## Secondary3800 0.74493 -0.15181 -0.59616 -0.9082469 -1.09817 0.45876
## Primary1800 1.24244 -0.54224 -0.06276 0.0002753 -1.28509 0.38255
## Primary2800 0.85424 0.18631 0.05840 0.3866843 0.60134 -0.32017
## Primary3800 0.55668 -0.75803 -0.15066 0.2617093 0.53754 -0.01956
##
##
## Biplot scores for constraining variables
##
## RDA1 RDA2 RDA3 RDA4 RDA5 PC1
## Slope 0.01146 -0.60054 -0.3615 0.04073 0.7120 0
## DistRiver -0.03892 -0.64987 -0.1478 0.04192 -0.7433 0
## soil_n -0.53844 -0.02599 0.1985 -0.75921 0.3060 0
## soil_p -0.54622 -0.02422 0.1950 -0.75312 0.3096 0
## soil_moist -0.16779 0.51504 -0.3419 -0.31827 -0.6989 0plot(trait.rda)